The nucleation, growth and stabilization of vortex lattices 



A. A. Penckwitt and R. J. Ballagh 

Department of Physics, University of Otago, Dunedin, New Zealand 

C. W. Gardiner 

School of Chemical and Physical Sciences, Victoria University, Wellington, New Zealand 

We give a simple unified theory of vortex nucleation and vortex lattice formation which is valid from the 
initiation process up to the final stabilization of the lattice. We treat the growth of vortex lattices from a rotating 
thermal cloud, and their production using a rotating trap. We find results consistent with previous work on 
the critical velocity or critical angular velocity for vortex formation, and predict the initial number of vortices 
expected before their self assembly into a lattice. We show that the thermal cloud plays a crucial role in the 
process of vortex lattice nucleation. 
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The opportunity that Bose Einstein condensates present for 
precisely controlled experimental and theoretical study of su- 
perfluids has revived an intense interest in the properties and 
behavior of quantized vortices. Many different theoretical cal- 
culations have been made of stationary single vortex states 
and their excitations ([[lj] and references therein), and vor- 
tex lattices (Q| and references therein) in Bose condensates, 
but the recent experimental observations of vortex lattices 
[^, BL [| ^, ^, [J have focussed attention on the mechanisms 
of vortex formation. 

One experimental method for creating vortices is to stir a 
condensate with an anisotropic potential, and a number of 
theoretical analyses of this scenario have been made in terms 
of the Gross-Pitaevskii equation, e.g. [ |To[ [T^ , 13 1. The 
view of these treatments is that the perturber causes a mixing 
of the condensate ground state and excited condensate states 
(with angular momentum values determined by the stirrer ge- 
ometry). Typically, a perturbative calculation is used to ob- 
tain a critical rotational speed (or a critical linear speed at 
the Thomas-Fermi radius) of the stirrer at which the mixing 
becomes effective. The argument is then made that an insta- 
bility will lead to growth of the vortex state. A more com- 
plete non-perturbative calculation for a localized rotating stir- 
rer [|l4|] shows that in this case the mixing is predominantly 
between the ground state and an I = 1 vortex state, and that as 
the relative amplitudes evolve by coherent "Rabi cycling", the 
full condensate exhibits a vortex cycling from infinity to the 
central regions of the condensate. No coherent mixing mech- 
anisms, however, can explain the formation of a vortex lattice, 
since an energy barrier exists between the superposition state 
and the vortex lattice state, and some additional mechanism is 
required to remove the energy liberated when the vortex lat- 
tice is formed. 

An alternative method of producing vortices, in which a va- 
por of cold atoms is evaporatively cooled so as to preserve 
its angular momentum has been demonstrated by Haljan et al. 
[H]. In this experiment, which involves no stirring potential, 
the mechanism of formation of the vortices must involve a 
modification of the theory of condensate growth (which now 
exists in a reasonably good quantitative form; [ |l5| ] and ref- 



erences therein) to take account of the non-zero angular mo- 
mentum of the vapor from which the condensate is formed. In 
this paper we will show that the mechanism thus demonstrated 
is the fundamental process in those experiments in which the 
stirring of the condensate apparently generates a vortex lat- 
tice. Weak stirring, in the absence of any thermal cloud, only 
produces "Rabi cycling" of angular momentum in and out of 
the condensate. However, if the stirring also produces a rotat- 
ing thermal cloud, then, by the same process responsible for 
the growth of a condensate [|T6||, angular momentum is trans- 
ferred irreversibly into the condensate, leading in equilibrium 
to a rotating vortex lattice. In general terms, this has already 
been noted in [|l7|], and in [|l8|] we derived a very simplified 
vortex growth equation — given below as Eq. (^f) — embodying 
the essentials of this concept. 

Other theoretical treatments of condensate dynamics at fi- 
nite temperature have been applied to the case of vortex dy- 
namics. Zhuravlev et al. [|l7| studied the interaction of vortex 
lattices with a thermal cloud in the presence of a trap asym- 
metry, and from broad conservation principles they obtained 
equations for the time evolution of the angular velocity of the 
condensate and thermal cloud, and the transfer of energy be- 
tween them. Williams et al. |19| showed how the inclusion 



of interactions with a thermal vapor cloud would give rise to 
an instability of surface modes, as required to produce vor- 
tices. Our own treatment Jl8| ] had a basis broadly equivalent 
to that of Williams et al., but showed how the interaction with 
the vapor cloud could produce an equation which should de- 
scribe the full process of vortex nucleation, growth and stabi- 
lization. Using this equation, we will show in this paper how 
the growth of vortex lattices can be understood as a conden- 
sate growth process, and that the condition that the gain in this 
growth process be positive is the same as the criteria noted 
above for vortex nucleation. We also demonstrate explicitly 
how the vortices produced in the "Rabi cycling" caused by a 
rotating trap nucleate the growth process, which then proceeds 
to dominate the vortex lattice formation process. 
The Gross-Pitaevskii equation for tpR = cxp(i fl • Lt/h)ip, 
the condensate wavefunction in a frame of reference rotating 
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with angular velocity ft, takes the form 

iK^=(-^V 2 + V* + u\Tp R \ 2 -n.L)Tp R , (1) 
at \ 2m ' ' J 

where V^?(x) is the trap potential in the rotating frame. A 
vortex lattice results from the stationary solutions of Eq. ([[]) 
satisfying ihdipu/dt = pcipR [[| which are normally 
obtained by integrating the Gross-Pitaevskii equation along an 
imaginary time direction [^0||. Indeed, noting that the achieve- 
ment of the vortex lattice stationary state requires dissipation 
in order to allow the system to settle into the energy minimum 
associated with the lattice, Tsubota et al. pl| ] proposed an ad 
hoc modification of this method as a phenomenological model 
of the process of vortex lattice formation. They modified the 
left-hand side of Eq. ([[]) to (1 + i^)ifidipn/ dt, and compen- 
sated for the consequent loss of probability by continuously 
adjusting the norm of the wavefunction. 
Origin of dissipation: When the condensate is grown from 
a rotating vapor cloud the growth process itself is dissipative, 
and in the experiments of [|J [|, ||, ||] dissipation can arise by 
transfer of atoms between the thermal cloud and the conden- 
sate. The vortex growth equation derived in [|l8| uses quan- 
tum kinetic theory for the growth of the wavefunction for a 
condensate trapped in a potential Vr(x, t) rotating with an- 
gular velocity ft, from a thermal cloud held at temperature T 
and chemical potential \i and rotating with angular velocity 
a. A dissipative term arises from collisions between atoms 
in a thermal cloud of atoms trapped by the same potential as 
the condensate, in which one of the colliding atoms enters the 
condensate after the collision. The net rate at which atoms en- 
ter and leave the condensate as the result of such collisions is 
characterized by a coefficient 

W+(x) = " / d 3 K 1( i 3 K 2 d 3 K3 5(a; 123 ) 
(27r) a 7i J 

xS{K 1 +K 2 -K 3 )F 1 F 2 (l + F 3 ) (2a) 
w g x 4:m(akT) 2 /nh 3 , (2b) 

and in practice g w 3 fits most growth experiments. Here 
F(x, K) is the Bose-Einstein distribution function for a ther- 
mal cloud in an effective potential in the frame rotating 
at the cloud's angular velocity a, i.e. V Q (x) = Vr(x) — 
m (a X x) /2, where V^r(x) is the time averaged trap po- 
tential over the differential rotation between the trap and the 
cloud. The vortex growth equation takes the form in the frame 
rotating with the trap potential 

ih ^tr = (-£v 2 + ^(x)+ u |^| 2 -o-l)^ 

+h yi> + (a - to) • L - iK^J i> R , (3) 

where 7 = KW + /kT. It is important to distinguish between 
the cloud's angular velocity a, and the trap 's angular velocity 
Q. If these are different the distribution of the thermal cloud 
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FIG. 1: Condensate density during formation of vortex lattice from 
a rotating thermal cloud. Vortices are marked by + or — according 
to their sense, (a) t = 29.28 to, (b) t = 37.70 to, (c) t = 52.40 to, 
(d) t = 187.24 to. Condensate wavefunction seeded as described in 
text. Parameters are a = 0.65 u), (i = 12hu>, 7 = 0.1, u = lOOOtto, 
initial nc = 12.77iw, Units are: time to = 1/u, distance ro = 
\/Ti/2muj, collisional strength uo = TuvrQ. 

is not truly thermodynamic equilibrium, and the form chosen 
represents an approximation expressing the experimental ob- 
servation that it is very difficult to spin up a thermal cloud 
by means of a rotating potential. One therefore expects that 
a < 51. In this paper, however, we shall consider mainly the 
cases in which either the trap is cylindrically symmetric, so 
that we may set ft = 0, or in which the cloud and the trap 
rotate at the same velocity, soq = tt. Eq. (Q) is essentially 
equivalent in terms of the physical assumptions required for 
its derivation to that of [|l9|], although its appearance is very 
different, and it is in practice easier to solve. In contrast, al- 
though Eq. ([3]) looks superficially very similar to that of [|2l|], 
it is in fact very different, as explained in Jl8| ] . 
Simulation results for a rotating vapor cloud: We have 
simulated Eq. (|3|) in two dimensions for many different sce- 
narios, and in Fig. [l] we present a time sequence of results 
from a representative case in which the condensate is ini- 
tially in the ground state of a rotationally symmetric harmonic 
trap of frequency u, the thermal cloud has p = 12huj, and 
a = 0.65 oj. The initial condensate has unit norm and chemi- 
cal potential nc = 12. 7 fad, although they can be any non-zero 
number. We simulate the non-stimulated collisions which start 
the process by adding to the intial wavefunction a uniform su- 
perposition of angular momentum states with I = 1 to 30, on 
a Gaussian radial profile centered at the Thomas-Fermi radius 
and with maximum amplitude ~ 2 x 10" 7 . 

Fig. |j](a) is a density plot of the condensate near the end of 
the first stage of the evolution, in which an imperfect ring of 
19 vortices arrives from infinity to just outside the Thomas- 
Fermi radius. As this ring shrinks further (Fig. |l](b)), several 
vortices are shed, and a dominant ring of 16 vortices passes 
through the Thomas-Fermi radius into the interior of the con- 
densate. These vortices then distribute themselves irregularly 
but quasi-uniformly over the condensate which expands and 
picks up angular velocity (Fig. |l](c)). Subsequently, over a 
long period, further vortices leave (and for larger values of a 
may enter) the dense region, until finally a regular lattice of 12 
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FIG. 2: Time evolution of condensate quantities for the case of 
Fig.[jJ (a) (L z ) and norm, (b) (/ii oc ) and <j m (in units of hui). 

vortices rotating at angular velocity a in the lab frame remains 
(Fig. id)). 

The overall time scale of the process is illustrated in 
Fig. ||(a) which shows the total angular momentum increasing 
sharply around t m 40, when the 16 vortices enter the dense 
part of the condensate. For a given seeding, the timescale 
for this onset is determined by 7 and to good approximation 
scales as 7 . The local chemical potential ^i oc (x) (defined 
as the absolute value of the right-hand side of Eq. (p divided 
by tpjf) provides a useful characterisation of the temporal de- 
velopment. In Fig. ||(b) , we plot the evolution of the spatial 
mean of /Ui oc ( x )> which exhibits a rapid initial adjustment (to 
w (j,), as the condensate roughly equilibrates with the thermal 
cloud, and a later dip at t w 48 during which the number of 
atoms in the condensate increases (see Fig. ||(a)). The spatial 
variance of jLti oc (x) quantifies the deviation of the solu- 
tion from equilibrium (where it is zero). It becomes large as 
the vortices enter the dense region of the condensate, and its 
subsequent slow decay indicates the long time needed for sta- 
bilization of the vortex lattice. It is not easy to characterize 
the latter time scale, since the local chemical potential can be 
almost constant everywhere, with possibly only a single vor- 
tex not in its final place. From a large number of simulations 
we have found that the number of vortices in the stable fi- 
nal lattice increases with a and /i, and is generally different 
from the number seen in the initial ring. It also exhibits a very 
weak dependence on the exact form of the seeding. For the 
case of fj, = YlTiuj the critical value for the appearance of any 
vortices is found from the simulations to be a = 0.444 ui, at 
which speed an irregular ring of 10 vortices appears initially, 
settling down eventually to three vortices in a lattice inside 
the Thomas-Fermi radius. Three dimensional simulations of 
Eq. (j|) also confirm the general behavior seen in Fig. |l| 
Analytic treatment: An analytic understanding of the nu- 
cleation and growth process of the vortices can be obtained 
from Eq. (|3|) by considering the linearized Bogolyubov exci- 
tations above the initial state, which have angular momentum 
I and energy eigenvalues e n ,i measured relative to nc- From 
the final term of Eq. (||) we see that excitations for which 
£n,z < M ~ He + hod will experience a positive growth rate 
G « j(fj, — fie + h&l — £n,l)/fi>- Initially this causes rapid 
growth of the condensate (the n = I = component, which 
has the lowest eigenvalue) until the condensate and thermal 
cloud approximately equilibrate (/i w fj,c)- Subsequently, the 
gain of the other components is j(al — e ni i/h), and thus the 
critical value of a for positive gain is given by eqj = hal, the 
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FIG. 3: (a) Evolution of angular momentum occupation probabili- 
ties Pi for I = 1 to 30 for the case of Fig. [lj (b) Comparison of gain 
coefficients of I components obtained from the simulation in Fig. |l| 
(solid line) with the predicted gain coefficients G = 7(0? — £q,i/K) 
(dotted line). 

same condition for criticality found in fll9|]. Above critical- 
ity, the dominant value l v is that which maximizes the gain, 
and is given by d{Kal — £o : i)/dl = 0, and a crude estimate 
for l v can be obtained by minimizing h 2 l 2 /2mr 2 — hal at 
the Thomas-Fermi surface, to give l v « nR 2 ?F 2ma/h. This 
is the number of vortices which would result from filling a 
disk of radius Rtf with a vortex lattice [E2|]. In the initial 
stages of vortex formation, when the process is still linear, the 
wavefunction of the condensate takes the form ( for simplicity 
including only the excitation with maximum gain, which has 
angular momentum Td v and energy Kui v = %al v — erj,Z„) 
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Here £o( r ) is the initial rotationally symmetric condensate 
wavefunction, <j> is the azimuthal angle, and u(r), v(r) are ob- 
tained by solving the Bogolyubov-de Gennes equations. The 
essential behavior can be seen by neglecting v(r), so that ip 
has l v zeroes — that is vortices — given by the interference of 
the two terms. These vortices all occur at the same radius, 
initially at infinity. Since the long distance behavior of u(r) 
must be less rapid than that of £,o(r), the ring will steadily 
shrink as the excitation grows. The time dependence of the 
coefficient means that the ring of vortices rotates at the an- 
gular frequency uj v in the rotating frame. In a more detailed 
picture, more values of I should be included, and the corre- 
sponding superposition will be an imperfectly circular ring of 
vortices, perhaps containing more than l v vortices. 

This analytic treatment can be verified by decomposing the 
condensate into angular momentum components, and plotting 
the occupation values Pi against time, as in Fig. [}](a). We see 
that the I components above threshold grow exponentially up 
to the time t « 40 (when the vortex ring passes through the 
Thomas-Fermi radius), and the gain is maximum for l v = 16. 
In Fig. ||(b) we compare the gain prediction j(al — eaj/fi) 
(obtained by calculating the surface modes eo,z for a conden- 
sate with fie = 12%u) ) against the gain rates measured from 
the simulation. Not only is our gain prediction very accu- 
rate, but we also find that the spatial particle density for each 
I component very accurately matches the prediction from the 
corresponding Bogolyubov wavefunction. 
Rotating trap, no thermal cloud: Finally, we consider the 
connection between stirring and growth. As an example of 
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FIG. 4: Evolution of angular momentum occupation probabilities 
Pi for: (a) rotating elliptical trap with no thermal cloud (7 = 0); (b) 
rotating elliptical trap and corotating thermal cloud. Only the even 
I states are shown. Parameters as in Fig. [j] except trap parameters 
uj x = 1.05 u>, (jj y = 1.15 to and a = fl = 0.65 u. 

a case of pure stirring we have simulated a rotating elliptical 
trap in the absence of a thermal cloud. The general behavior is 
well represented by the case fie = 12.7fiw and il = 0.65 w, 
for which we find that the condensate deforms into an ellipti- 
cal shape rotating at the trap frequency, and a similarly de- 
formed ring of about 28 vortices shrinks until it meets the 
dense part of the condensate at its narrowest radius. The re- 
sulting state is a distorted version of Fig. 0(b) [|Z3||, and this 
subsequently undergoes an oscillatory behavior, but with no 
further penetration of the vortices into the condensate. Only 
the even angular momentum components of the condensate 
become significantly occupied and undergo periodic cycling 
(Fig. Q(a)). This can be well understood in terms of the Rabi 
cycling model fll-i]], noting that an elliptical potential can con- 
nect the initial ground state only to states of even I. Resonant 
mixing to the surface modes will occur when Mil — eo,i = 0, 
and the lowest value of ft for which this occurs (when also 
d(fiQl — £Q t i)/dl = 0) gives the same condition for criticality 
as found in Refs. [Ill, 12 13|], This is also exactly the same 
condition as for the rotating cloud. 

Rotating trap with a thermal cloud: Now repeating this 
simulation, but with a thermal cloud (Eq. (|^) with a = SI 
and 7 7^ 0), we get dramatically different results. The same 
elliptical distortion of the condensate occurs, but the effect 
tends to smooth out with time. A distorted vortex ring ap- 
pears, with an initial number of vortices as expected from the 
growth mechanism (here w 14), and then shrinks into the con- 
densate and nucleates the vortex lattice. The angular momen- 
tum projections in Fig. ^b) show the difference very clearly. 
The vortex lattice is the result of the growth process, which 
is essentially superposed on the effects of stirring, eventually 
masking them. This scenario leads to the conclusion that the 
principal effect of the rotating trap is to produce a rotating 
thermal cloud, which then transfers angular momentum into 
the condensate, eventually forming the vortex lattice. 
Comparison with other work: The main point of this pa- 
per is the analysis of vortex formation in terms of the growth 
in occupation of modes with non-zero angular momentum. 
This requires the presence of a rotating thermal cloud, but stir- 
ring the condensate can seed angular momentum components 
which may then grow from the thermal cloud by stimulated 
collisions. The critical angular velocity of the thermal cloud 
to allow growth is the same as the critical angular velocity re- 



quired for stirring to resonantly mix surface modes with the 
condensate, and agrees with the latter values calculated by 
previous workers g However, our treatment 

goes beyond the critical point, after which the behavior of the 
two processes becomes quite different, and includes the full 
dynamics of lattice stabilization. Our vortex growth equation 
(|3|) can be solved without linearization or using hydrodynamic 
methods, and its solutions can be interpreted in a straightfor- 
ward quantum mechanical way. The dynamics of the thermal 
cloud are suppressed in this paper, but their inclusion is nec- 
essary for the full understanding of vortex lattice dynamics, 
including lattice decay, which can also be treated by our for- 
malism. 
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